home *** CD-ROM | disk | FTP | other *** search
/ El Mac 9 / El Mac 9.iso / Shareware / Demos / Igor Demo Pro / 1 PutContentsIn Igor Pro Folder / WaveMetrics Procedures / Analysis / DSP (Fourier etc) / FTMagPhaseThreshold < prev    next >
Encoding:
Text File  |  1994-06-03  |  2.8 KB  |  94 lines  |  [TEXT/IGR0]
  1. #include <BringDestToFront>
  2. #include <Math Utility Functions>
  3.  
  4. | Given an input wave ( doesn't have to be a power of two), create a new wave
  5. | with the suffix "_Mag" that contains the normalized frequency response and
  6. | optionally a wave with the suffix "_Phase".
  7. | Several levels of resolution enhancement (really just sin x/x interpolation) are provided.
  8. | Options include windowing (Hann) vs no windowing, linear vs dB, phase vs no phase
  9. | If phase then option of radians or degrees,wrapped or unwrapped
  10. | The thresholdPct parameter suppresses phase
  11. | values where the magnitude is less than some percentage of the maximum magnitude.
  12. |
  13. | You may want to modify the code at the end of this macro.  It sets the display
  14. | style and that is a matter of taste.
  15. |
  16. Macro FFTMagPhaseThreshold(w,window,resolution,linlog,phase,phasetype,thresholdPct)
  17.     string w
  18.     Prompt w,"Input data:",popup WaveList("*",";","")
  19.     variable window=1
  20.     Prompt window,"Windowing:",popup "None;Hann"
  21.     variable resolution=1
  22.     Prompt resolution,"Resolution enhancement:",popup "none;2;4;8;16;32"
  23.     variable linlog= 2
  24.     Prompt linlog,"Magnitude Mode:",popup "Linear;dB"
  25.     Variable phase= 1
  26.     Prompt phase,"Phase:",popup "No phase;Phase in radians;Phase in degrees"
  27.     Variable phasetype=1
  28.     Prompt phasetype,"Unwrap phase?",popup,"No;Yes"
  29.     Variable thresholdPct=5
  30.     Prompt thresholdPct,"Suppress phase if below % max magnitude"
  31. ;
  32.     PauseUpdate; silent 1
  33.     
  34.     string destw=w+"_Mag",phasew= w+"_Phase"
  35.     Variable n= numpnts($w)
  36.     
  37.     if( (resolution<1) %| (resolution>6) )
  38.         Abort "resolution out of range"
  39.     endif
  40.     resolution -= 1
  41.     Duplicate/O $w $destw
  42.     if(window==2)
  43.         Hanning $destw; $destw *= 2            | assumes continuous rather than pulsed data
  44.     endif
  45.     Redimension/N=(CeilPwr2(n)*2^resolution) $destw        | pad with zeros to power of 2
  46.     fft $destw
  47.     $destw= r2polar($destw)
  48.     | NOTE: depending on your application you may want to un-comment the next line
  49. |    $destw[0] /= 2                                | dc is special
  50.     if( phase!=1 )
  51.         Duplicate/O $destw $phasew
  52.         Redimension/R $phasew
  53.         $phasew= imag($destw)
  54.         | remove phase values where corresponding amplitude is small    | FLORIN
  55.         if( thresholdPct > 0 )
  56.             WaveStats/Q $destw
  57.             Variable threshold = V_max * thresholdPct/100
  58.             $phasew*= ($destw[p] > threshold )
  59.         endif
  60.         if( phasetype==2 )
  61.             $phasew[0]= $phasew[1]            | try to avoid glitch at dc
  62.             UnWrap 2*Pi,$phasew
  63.             $phasew[0]= 0
  64.         endif
  65.         if(phase==3)
  66.             $phasew *= 180/Pi
  67.             SetScale y,0,0,"deg",$phasew
  68.         else
  69.             SetScale y,0,0,"rad",$phasew
  70.         endif
  71.     endif
  72.     Redimension/R $destw
  73.     if( linlog==2 )
  74.         WaveStats/Q $destw
  75.         $destw= 20*log($destw/V_max)
  76.         SetScale y,0,0,"dB",$destw
  77.     else
  78.         $destw /= n/2
  79.         SetScale y,0,0,"V",$destw
  80.     endif
  81.     BringDestFront(destw)
  82.     if( phase!=1 )
  83.         CheckDisplayed  $phasew
  84.         if( !V_Flag )
  85.             Append/R $phasew
  86.         endif
  87.     endif
  88.     if( numpnts($destw) <= 129 )
  89.         Modify mode($destw)=4,marker($destw)=19,msize($destw)=1
  90.     else
  91.         Modify mode($destw)=0
  92.     endif
  93. end
  94.